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Abstract 



A recent theorem giving the initial behavior of very short-time fluctuations of 
particle displacements in classical many-body systems is discussed. It has applica- 
{Sj ■ tions to equilibrium and non-equilibrium systems, one of which is a series expansion 

ON , of the distribution of work fluctuations around a Gaussian function. To determine 

the time-scale at which this series expansion is valid, we present preliminary numeri- 
cal results for a Lennard-Jones fluid. These results suggest that the series expansion 
converges up to time scales on the order of a picosecond, below which a simple Gaus- 
sian function for the distribution of the displacements can be used. 

Keywords: Fluctuations, cumulants, work, displacements, time expansion, non-Gaussian 
effects, Van Hove self-correlation function, Green's functions, Lennard-Jones fluid. 
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1 Introduction 

> ■ 

^ . In recent years, non-equilibrium fluctuations in entropy production, work and heat have 

become the subject of intense study. In this context, results such as the stationary state 
fluctuation theorem [1] have given a deeper insight in the behavior and character of non- 
equilibrium systems beyond irreversible thermodynamics and hydrodynamics, which are 
especially important for small systems such as biomolecules and nanotechnology. 

The fluctuation theorem of a non-equilibrium steady state considers fluctuations at 
long time-scales. Fluctuations of work, heat and entropy production at short time-scales 
are however less often studied. Such fluctuations play a role in short-time (picosecond) 
non-equilibrium transport at the nanometer scale such as in incoherent neutron scattering, 
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among others. While fluctuations of many properties can be studied, the present paper will 
give an overview of recent results for short-time fluctuations of the displacements of indi- 
vidual particles in a classical-mechanical system. In particular, we will focus on a theorem 
regarding the probability distribution of the short-time fluctuations of displacements [2j [3]. 
As we will show in this paper, this distribution is also related to that of work fluctuations. 
Furthermore, we will discuss its implications and range of applicability. 

The paper is structured as follows: In section [2] we introduce moments and cumulants 
of a probability distribution, in terms of which the theorem is formulated. In section [3] 
we give a simple example of the short-time behavior of cumulants as a prelude to the 
theorem presented in section HI One equilibrium and one non-equilibrium application of 
this theorem are discussed in section [51 Preliminary results on the range of validity of the 
short-time expansion and thus of the utility of the theorem are presented in section [61 We 
end with the conclusions in section [71 



2 Moments and cumulants 

Consider a classical system of N particles with positions r$ and velocities Uj. The equations 
of motions are Hamiltonian with forces derived from a smooth potential U(r\, . . . rjv), while 
the initial phase points have a probability distribution denoted by p{{ri, In the current 
context, Lennard- Jones potentials may be considered smooth, despite the singularity at 
zero, provided initial conditions which hit this singularity have measure zero. 

In time, the particles will be displaced from their initial position by an amount Arj(t) = 
Tiit) — Tj(0) = (Axi(t), Ayi(t), Azi(t)). The probability distribution of these displacements, 
P(Ari, . . . Avn, t), is determined by the dynamics and the initial probability distribution. 

Instead of considering the probability distribution function P(Ari, . . . Ar^, t) itself, it 
is often useful to consider the moments of this distribution which are defined as 

H ni ... nN (t) = (Ax^Ay^AzT* • • • As^-AyJ^A^')*, (1) 

where the average ( ) t is taken with P(At*i, . . . Ar^, t) and each rij is triplet of integers, 
rti = (ni X ,niy,ni Z ). These moments determine the moment generating function 

p( kl ,...k N ,t)- i^—j- \-^r \~^T'" \-^r^- nAt)t 

ni x =0 ni y =U y ni z =U "jv z =U 

(2) 

which coincides with the Fourier transform of P(Ar*i, . . . Arjy, t). Thus all information 
about P(Ari, . . . Ar^r, t) is in principle contained in these moments. 

A special case is formed by taking rii = (n, 0, 0) and all other n^i = 0, which means 
one considers the distribution of the displacement of a single particle, P(Ax,t). Its mo- 
ments are given by /x„(t) = (Ax") t and its generating function by P(k,t) = J^^Lo Vnif) — 
(exp(ikAx)) t - In equilibrium, P(k,t) becomes equal to the incoherent intermediate scat- 
tering function that can be measured e.g. by neutron scattering. 
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In special cases, such as for an ideal gas or a harmonic lattice in equilibrium, the moment 
generating function is Gaussian and then determined by the first and second moments, i.e, 
for the single particle displacement by the p\(t) and p,i{t). Even in those cases, the higher 
order moments are not zero, but are simple factored forms involving \i\ and p 2 - 

A more convenient representation of the distribution for Gaussian and near-Gaussian 
distributions is in terms of cumulants, which are, loosely stated, moments with all possible 
factorizations taken out. Formally, they are defined as the derivatives of the logarithm of 
P, such that: 

{iky 



n=l 



-K n (t) =logP(M). (3) 



This equation defines the cumulants K n (t) for the single particle displacement distribution 
P(Ax, t), but it can easily be generalized to deal with the general case of P(Ar 1 , . . . Ar N , t). 
Note that in case of a Gaussian distribution, the cumulants K n with n > 2 are all equal to 
zero. 

One can find explicit expressions for the cumulants by substituting the expression on 
the right-hand side of equation <^ into the right hand side of equation expanding the 
logarithm and equating terms of equal powers in k. For the first four powers of k this 
yields 

Ki{t) = (4a) 
K 2 {t) = p 2 (t)-p\{t) (4b) 
k 3 (*) = // 3 (*)-3pi(*)/z 2 (*) + 2nl(t) (4c) 
iH(t) = /i4(t)-4^i(*)/i 3 (*)-3^(t) + 12^(t)/x 2 (*)-6/x*(*). (4d) 

These first four cumulants are also known as the mean, variance, skewness and kurtosis, 
respectively. 



3 Short-time behavior: a harmonic example 

As an illustrative example, consider a single particle of mass one in one dimension, subject 
to a harmonic potential U = x 2 /2. Let the initial velocity distribution be Gaussian: 

p{v ) = (2vrA; B T) 1 / 2 exp 

The distribution of the initial positions xo will be arbitrary and given by some function 
p(xq). Because p is arbitrary, this single-particle system is not necessarily in equilibrium. 
The equations of motion of this system are easily solved and give 

x(t) = cost xq + sin tvo. (6) 



2k B T 



(5) 
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Figure 1: Short-time behavior of cumulants of the displacement of a particle in a harmonic 
potential, initially not in equilibrium, based on equation (8). The left plot shows the 
cumulants function of time t, for n = 2, 3, 4, 5 and 6, respectively. The right plot 

shows the same on a log-log scale, to bring out the scaling behavior for small t. 

Using equations fl^D and (jSJ), the function P(k,t), i.e., the Fourier transform of P(Ax,t) 
with Ax = x(t) — x , can now be computed: 



P(k,t) = ldx jdv p(x )p(v y k ^-^ 

dx J dv p(xo)p(v )e ik[{cost - 1)xo+s ' mtvo] 
dx p(x )e iHcost - l)xo [ dv p(v )e lksintvo 
p(k[cost — 1]) exp 



— KflTsin t 
2 



(7) 



where p is the Fourier transform of the distribution of initial positions, p(xq). Given this 
result for P(k,t) and equation ([3]), one finds for the cumulants: 

Kn (t) = 4 0) (cost - l) n + 5 nt2 k B T sin 2 1, (8) 

where Kn^ is the nth cumulant of p(x ). 

An example of the cumulant behavior in equation (jSj) is given in figure [TJ For this 
example, we chose k B T = 1 and for p(x ) a non-Gaussian function which is continuous and 
piece-wise linear: 

{|(l + x ) if - 1 < x < | 
|(1 -x ) if | < x < 1 
otherwise. 

The reason for choosing this p(xo) is that its cumulants are non-trivial yet easy to compute. 
The first panel in figure [1] shows the behavior of the cumulants n n for n = 2, 3, 4, 5 and 6 
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as a function of time. The second panel shows the same quantities, but zooms in on the 
small time behavior and has a log-log scale. On this log-log scale, the curves for n = 2, 3, 
4, 5 and 6 have, for small t, slopes —2, —6, —8, —10 and —12, respectively. This means 
that for small t, «2 oc t 2 , while «2<n<6 oc t 2n . It can easily be verified that equation (jHJ) 
supports this observation for all choices of p(xq) and for all n > 2. 

One may wonder whether the t 2n behavior of K n (t) for small times t is general, or 
whether it is due to the simplicity of our example. In fact, this behavior is more general, 
under conditions explained in the next section. 

4 Theorem on short-time fluctuations of displacements 

Cumulants of particle displacements have been studied in the context of the intermediate 
scattering function in the field of neutron scattering on equilibrium systems for a long 
timejH, El IB]- It was found that while all odd cumulants of displacements are zero in 
equilibrium, the even cumulants for small times obey K n (t) oc |t| n+1 for hard spheresjB], 
while for systems with smooth potentials it was known that K±(t) oc t 8 [I] and Ke(t) oc t 12 



Unfortunately, these results for smooth potentials only held for equilibrium systems and 
furthermore no result s for general n were available. A recent theorem however, generalizes 
these results to certain non-equilibrium systems and to general n[21 E]- It states that if 

1. the forces are smooth functions of the positions of the particles (smooth here means 
that all derivatives exist for all points in phase space except a set of measure zero), 

2. the forces are independent of the velocities, 

3. the initial velocities are Gaussian distributed, and 

4. the initial velocities are independent of the initial positions, 
then the cumulants of particle displacements satisfy for small times t 



where the c n are constants. In fact, all n-order cumulants satisfy this rule (an n-order 
cumulant is cumulant K ni nN such that ^(n^ + n iy + n iz ) = n). 

We stress that this theorem applies to non-equilibrium as well as to equilibrium systems, 
because the (initial) position distribution can be arbitrary and does not have to be equal to 
its equilibrium form. A common class of non-equilibrium systems covered by the theorem 
consists of those that start out in equilibrium and are then driven out of equilibrium by 
some external perturbation, because the requirement of a Gaussian velocity distribution 
only has to be satisfied at the initial time. 

A sketch of the proof of this theorem was given in reference 2, while the full proof can 
be found in reference 3; we will therefore not prove the theorem here again. However, we 



15]. 




(9) 
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do want to remark that the nature of the dynamics is largely unimportant for the theorem 
to hold (as long as the forces are smooth), and that the most important requirement seems 
to be that the velocities are initially Gaussian distributed. To what extent one can relax 
the requirement that this Gaussian distribution be independent of the initial positions, is 
at present not clear. 



5 Applications 



5.1 Neutron scattering in equilibrium systems 

As mentioned above, the short time behavior of the cumulants of displacement arises in 
the context of neutron scattering. In fact, in equilibrium, the incoherent structure function 
F s (k, t) is precisely the Fourier transform of the distribution of single particle displacements 
P(k,t)^\. Given equation Q, we can therefore write 



F s (k,t) = exp 



(10) 



Note that all odd cumulants are zero for an equilibrium system for symmetry reasons. 

For small wave- vector k, the terms involving higher cumulants in equation fllOp can 
be neglected and F s (k,t) becomes approximately Gaussian, as does its Fourier inverse, 
the Van Hove self-correlation function G s (Ax,t). While equation (flOl) in principle allows 
for a systematic improvement on this Gaussian approximation by increasing the number 
of terms considered in the cumulant sum, in its current form, it is not possible to take 
into account higher order cumulants and still perform the Fourier inverse analytically. To 
obtain a systematic analytic improvement upon the Gaussian approximation of the Van 
Hove self-correlation function one first has to expand the right-hand side of equation (flQj) 
around a Gaussian as follows: 



F s (k,t) = exp 



~^2(t)k 2 



1 + 



«4(0 ,4 _ K eW 



4! 



6! 



k 6 + 



11^ 



By taking the Fourier inverse of this expression, we find for the Van Hove self-correlation 
function 



G s (Ax,t) 



exp(- 



-w 



"2-KK 2 {t) 



K 4 (t)H 4 (w) K 6 (t)H 6 (w) 
4!4k§(*) 6!8«ij(*) 



+ 



'121 



where w = Ax / \/2k 2 and H n is the nth Hermite polynomial. 

Provided that the system has smooth forces independent of the velocities (conditions 
1 and 2 of the theorem, respectively), the third and fourth conditions are automatically 
satisfied in (canonical) equilibrium. The application of the theorem now consists of noting 
that, with w = 0(1), the terms inside the brackets on the right-hand side of equation ffl2|) 
are 1, t 4 , t 6 , respectively. It can be shown that this behavior persists, i.e., each next term 
is two powers of t higher than the previous one. Hence as a consequence of the theorem, 
the subsequent terms in this expansion become successively smaller if t is small, so that a 
truncation of the expansion can give an accurate approximation to G s (Ax,t). 
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5.2 Distribution of work done by an external field 



Since the theorem is not restricted to equilibrium systems, we want to discuss at least 
one interesting application of the theorem to a non-equilibrium system, namely a system 
in a constant external field along the x direction. This constant external field exerts a 
constant force F on all N particles, which does an amount of work W on the system. Since 
the initial conditions are drawn from a probability distribution, W will be a fluctuating 
quantity. Non-equilibrium fluctuations of work have been investigated in the context of the 
fluctuation theoremp] and the non-equilibrium work relation[8]. The theorem of section H] 
can be used to investigate a different aspect of the work fluctuations, namely their very 
short-time behavior. 

Consider now the work W done on the system by this force during a time t: 



W = Y] f dt' F-Viit') 
i=i Jo 

N 

= \F\^Axi(t). 



i=l 



Given this linear relation between the displacements and the work, the short-time theorem 
can be applied to the distribution of work values as well. To see this, one can apply the 
multinomial expansion[3] to the nth cumulants of the work distribution, i.e. 



K n(t)~ E ni Wsn N \ Knin2 * nN{t) ' (13) 

where n nin2 s nN = K nin2 TlN with rij = (n,-, 0, 0). Examples are = iV«ioo... an d ^ = 
NK200... + N(N — l)«noo... (where we used that the particle are indistinguishable). If the 
forces are smooth and independent of the particles' velocities and if the initial velocities 
of the particles are Gaussian distributed, then according to the theorem, each term in the 
sum on the right hand side of equation (fT3"|) is of order t 2n if n > 2. In other words, the nth 
cumulant of the work fluctuations (t) also scales as t 2n for n > 2 if t is small enough. 
Consequently, the distribution of work fluctuations may be written in a similar expansion 
as that of the Van Hove self-correlation function in equation (112p . 



6 Relevant time scale for an equilibrium Lennard- 
Jones fluid 

We have discussed the short-time behavior of the cumulants of displacements, but we 
have not yet addressed the question of how short the time scale has to be for the above 
applications to be meaningful. We will now address this question now for single particle 
displacements, which play a significant role in the Van Hove self-correlation function G s . 
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We would call the series expansion of G s in equation (fl2|) meaningful if each next term is 
smaller than the previous one. While the theorem guarantees that this is so below some 
time-scale, it does not say what this time scale is. An estimate for the time scale at which 
these short-time results will no longer be applicable could be found by considering the 
time t = tq at which the first correction term in equation ( fl2l) becomes of comparable 
magnitude as the first term (the Gaussian). In that case, the series in equation ffl2|) can 
no longer be expected to be useful. Note that in equation ffT2l all odd cumulants have 
vanished because of time-reversal symmetry, so the first correction term involves k 2 and 
^4, or, using equation ([9]), the values of the constants c 2 and C4. Since w = 0(1), one can 
find tg by equating the first and (the absolute value of) the second term in equation ( [121) 
with t set to r G and H{w) set to 1: 



Expression for c n can be found in section 3.4 of reference 3. While the values of the c n 
will be system dependent, but it is nonetheless of interest to have a realistic estimate of tq- 
For this purpose, we performed molecular dynamics (MD) simulations at constant volume 
and energy (NVE) for an equilibrium iV = 100 particle system with periodic boundary 
conditions |9j. The inter-atomic potential used is the Lennard- Jones (LJ) potential V{r) = 
4e[(a/r) 12 — (cr/r) 6 ]. All quantities reported are in Lennard- Jones units: length in units of 
a, temperature in units of e/Icb, number density (p) in units of cr -3 and time in units of 
= (cr 2 m/e) 1 / 2 , where m is the mass of the particle. Since these are arbitrary units, to 
understand the physical consequences of our results, we use the LJ parameters of a specific 
substance as a reference, namely Argon. In that case, t lj corresponds to 2.16 picoseconds 
while e/ks = 119.8 K [3]. In the simulation, a potential cutoff of r c = 2.5a was used 
and the equation s of motion were integrated using the Verlet algorithmic] with a time 
step of one femtosecond. Data were accumulated once equilibrium had been attained in 
the simulation and collected between sufficiently long periodic time intervals, in order to 
ensure statistically uncorrelated data points. 

Note that while strictly the theorem should be applied to a system in the canoni- 
cal (NVT) ensemble (to ensure a Gaussian velocity distribution), the preliminary data 
presented here are obtained in the microcanonical ensemble (NVE). Since in the thermo- 
dynamic limit, these two ensemble should become equivalent, the time scale tg obtained 
in the NVE simulation is expected to be very close to the tq that would be found in an 
NVT simulation. 

In figure (2J we show preliminary results from these simulations for C4 as a function of 
the kinetic temperature T. More results will appear in reference 10. The value of c 2 = ksT, 
so from the results in figure [2] we can estimate the time scale tq using equation ffl4")) . It 
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c 4 t g 



4!4[c 2 t 2 P 



where equation ([9]) was used. Solving this equation for tq yields 




(14) 
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Figure 2: Numerical results for C4, the coefficient in front of the t 8 behavior of ^4, as a 
function of the kinetic temperature T. The horizontal lines through the points are the 
error-bars. The system is a Lennard- Jones fluid with density p = 0.8. These results are 
from a MD simulation in an NVE ensemble with N = 100 particles, with periodic boundary 
conditions. All quantities are in LJ units. 

is found that this time scale ranges from tq = 0.3 for T = 1 to tq = 3 for T = 3, all 
in Lennard- Jones units. This corresponds to time scales on the order of 0.6 picoseconds 
for T = 119.8 K to 6 picoseconds for T = 359.4 K. We note that these time scales are 
to be interpreted as upper bounds on the times at which expansions such as those in 
equation ffl2l break down. 

7 Conclusions 

We discussed a recent theorem which states that for smooth systems with initial Gaussian 
velocity distributions (but with arbitrary position distributions), the cumulants satisfy a 
scaling relation for short times: K n oc t n if n < 3 and K n oc t 2n for n > 3. As applications 
of this theorem, this scaling implies that the expansion of the Van Hove self-correlation 
function around a Gaussian is useful, as is a similar expansion for the short-time work 
fluctuations in a non-equilibrium system subject to a constant external field (provided the 
initial distribution of velocities is Gaussian). 

We also showed some preliminary numerical results for a Lennard- Jones model of Argon 
which suggest that this expansion breaks down for larger times, with a n upper limit of 0.6 
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to 6 ps depending on the kinetic temperature. Further details will be published in a later 
paper [TO]. 

Finally, we wish to note that other applications of the theorem are the scaling of non- 
Gaussian parameters characterizing "dynamical heterogeneities" in supercooled liquids and 
glasses, and non-equilibrium transport using a Green's functions approach; see reference 3 
for details. 
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